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ABSTRACT 

We demonstrate that the output of a cosmological N-body simulation can, to remark- 
able accuracy, be scaled to represent the growth of large-scale structure in a cosmology 
with parameters similar to but different from those originally assumed. Our algorithm 
involves three steps: a reassignment of length, mass and velocity units, a relabelling 
of the time axis, and a rescaling of the amplitudes of individual large-scale fluctu- 
ation modes. We test it using two matched pairs of simulations. Within each pair, 
one simulation assumes parameters consistent with analyses of the first-year WMAP 
data. The other has lower matter and baryon densities and a 15% lower fluctuation 
amplitude, consistent with analyses of the three-year WMAP data. The pairs differ 
by a factor of a thousand in mass resolution, enabling performance tests on both 
linear and nonlinear scales. Our scaling reproduces the mass power spectra of the 
target cosmology to better than 0.5% on large scales (fc < Cl.lft.Mpc~ 1 ) both in real 
and in redshift space. In particular, the BAO features of the original cosmology are 
removed and are correctly replaced by those of the target cosmology. Errors are still 
below 3% for fc < 1/iMpc -1 . Power spectra of the dark halo distribution are even 
more precisely reproduced, with errors below 1% on all scales tested. A halo-by-halo 
comparison shows that centre-of-mass positions and velocities are reproduced to bet- 
ter than 90ft~ 1 kpc and 5%, respectively. Halo masses, concentrations and spins are 
also reproduced at about the 10% level, although with small biases. Halo assembly 
histories are accurately reproduced, leading to central galaxy magnitudes with errors 
of about 0.25 magnitudes and a bias of about 0.13 magnitudes for a representative 
semi-analytic model. This algorithm will enable a systematic exploration of the cou- 
pling between cosmological parameter estimates and uncertainties in galaxy formation 
in future large-scale structure surveys. 

Key words: cosmology:theory - large-scale structure of Universe. 



1 INTRODUCTION 

Ever since the first measurements, made almost 40 years 
ago, the clustering of galaxies has been recognized as a 
powerful and robust tool to investigate some of the most 
fundamental questions in cosmology and galaxy formation. 
Subsequent years have produced major advances on both 
observational and theoretical fronts, and the study of the 
three-dimensional distribution of galaxies has contributed 
enormously to our current understanding of cosmic evolu- 
tion. For instance, it has provided information about the to- 
tal matter content of the Universe, the amount of baryons, 
the formation and evolution of galaxies, and even about the 
nature of gravity. 

The potential to learn from the galaxy distribu- 
tion is still far from exhausted. New galaxy surveys 
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are scanning ever larger volumes providing more homo- 
geneous and more detailed structural information, and 
thereby improving considerably the accuracy which can 
be achieved with existing probes of cosmological physics. 
In addition, exciting new probes such as the imprint of 



Baryonic Acoustic Oscillations (BAO, e.g. [E iscnstci n et al 
20051 ; ICole et al.|[2005l ; iGaztanaga et al] 2008 ; Sa nchez et al 



2009), the form of redshift-sp ace distortions (jGuzzo et al 
20081 - IPercival fc White! 120081). the c l ustering amplitude 



of d a rk matter haloes l|Seliakl 120081 ; IMatarrese fc Verdd 
2008; ICarb one et all l2008h and the evolu t ion of gravita- 
tional lensing signals ( Huterer et al.l 200(: j iHoekstra et al.l 



2006; iBerge et alj|2009l ; iDebono et al.l 



2009) have been pro- 



posed to explore, for instance, the nature of the dark en- 
ergy and the amplitude of primordial non-Gaussianities. 
In addition, higher-order clustering could provide fresh in- 
sights into many of the cosmological issues mentioned ab ove 
l|Baugh et alJl2004lCr"oton et al.ll2007l ; iRoss et al.ll2008T l. 
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A key step towards achieving optimal exploitation of 
galaxy clustering data and towards assessing the power of 
future experiments is to refine our understanding of the con- 
nection between galaxy formation and the background cos- 
mological model. In particular, there are two aspects that 
need to be addressed. The first is the interplay between local 
baryonic physics an d the large-scale evolution o f density per- 
turba tions (see e.g. IZheng fc Weinberg! 120071 ; lHarker et al.l 
120071 ). For instance, how different is the efficiency of star 
formation in a universe with larger dark matter perturba- 
tion amplitudes? Can we use the small-scale clustering of 
galaxies, which is heavily affected by astrophysics, to learn 
about the underlying cosmological model? Another crucial 
aspect is to study whether uncertainties in galaxy formation 
physics will weaken or possibly bias the cosmological con- 
straints extracted from future datasets, e.g. how robust is a 
measurement of the dark energy equation of state derived 
from BAO to uncertaint ies in how (and where) galaxies form 
(e.g. lAngulo et aill2008l )? 

In order to explore these issues we require accurate 
and flexible theoretical predictions for the properties and 
the spatial distribution of galaxies over large volumes. Cur- 
rently, the most sophisticated models are built in a three- 
step process. First, the abundance, structure and cluster- 
ing of haloes are predicted by following the evolution of 
the dark matter in a h igh-resolution N-body simulation 
fe.g. lSpringel et ahlbood ). Second, the properties of galaxies 
within these haloes are p redicted using either semi-analytic 
recipes (e.g. Baugh 2006) or empirically motivated paramet- 
ric techniqu es known as halo occup ation distribution (HOD) 
models (e.g. lCoorav fc Shetlj|2002]) . Lastly, the appropriate 
flux limits, sample selection criteria, redshift completeness 
and survey geometry are applied to build realistic mock cat- 
alogues. 

Unfortunately, the realism of these mock catalogue im- 
plies a substantial computational cost, especially when very 
large scales are probed. A simulation that resolves simulta- 
neously Gigaparsec scales and Milky- Way sized haloes takes 
at least several days in the largest supercomputers avail- 
able today. The modelling becomes even more expensive 
when one considers that a large number of independent N- 
body simulations are needed to assess clustering errors prop- 
erly. An adequate estimation of the variance requires sev- 
eral dozen realizations of the density field (e.g. a 10% error 
on a variance estimate for a Gaussian distribution requires 
~ 50 realizations), and an order of magnitude more simu- 
lations are needed to robustly compute the full covariance 
matrix (jTakahashi et al]|2009l ). In practice, these technical 
constraints mean that this method can be applied to simu- 
late galaxies only in a very limited number of cosmological 
scenarios, which then very sparsely sample the allowed cos- 
mological parameter space. This results in an incomplete 
picture of the coupling between galaxy formation and cos- 
mology. 

In this paper we propose and test a new approach 
which avoids this limitation and allows high-fidelity mock 
catalogues to be created for a broad range of cosmological 
scenarios. Current constraints on cosmological parameters, 
derived from analyses of CMB fluctuations, of supernovae 
and of galaxy clustering, allow only relatively small varia- 
tio ns around the prefer red structure formation model (see 
e.g ISanchez et all I2009T ) . Taking advantage of this, we ex- 



plore whether it is possible to manipulate the results of a 
single N-body simulation to represent the growth of struc- 
ture (with the accuracy needed for faithful galaxy modelling) 
over the whole allowed cosmological parameter space. In par- 
ticular, we will focus on reproducing the position, velocity, 
abundance and mass of the dark matter haloes thought to 
host galaxies. Matching these properties accurately at all 
times ensures a correct merger history for all DM haloes 
and therefore the correct prediction of the properties of the 
galaxies within them for any given galaxy formation model. 
The resulting "scaled" galaxy catalogues thus correspond to 
those expected in the target cosmology. 

Our proposed algorithm consists of the following steps. 
First, a scaling of the box-size and particle mass together 
with a reassignment of the output redshifts of the original 
N-body simulation matches the shape of the linear variance 
of fluctuations in the target cosmology over the range of 
scales which give rise to observable nonlinear structure in 
that cosmology. Second, using linear perturbation theory, 
the large-scale components of the displacement and veloc- 
ity fields are corrected to take account of the differences in 
power spectrum shape between the original and target cos- 
mologies. 

The layout of our paper is as follows. In §2 we present 
more details of our algorithm. We discuss how length, mass 
and velocity units are altered and the time axis relabelled 
to match the nonlinear evolution of perturbations. We then 
describe how to rescale the amplitude of individual long- 
wavelength modes at each redshift in order to correct the 
linear and quasi-linear evolution on these scales. In §3.1 we 
present a suite of N-body simulations specially designed to 
test the accuracy of our scaling procedure and in §3.2 we 
use them to show how well the full nonlinear matter power 
spectrum and its evolution can be reproduced for a chosen 
cosmology by scaling the particle distribution from a simu- 
lation of a different cosmology. The resulting halo catalogues 
are then compared in §4 with catalogues constructed from 
a matched simulation carried out directly in the target cos- 
mology. Finally, in §5 we provide a summary and we discuss 
applications of our method. 



2 SCALING AN N-BODY SIMULATION 

In this section we will show how it is possible to use the re- 
sults of an N-body simulation carried out in one cosmology 
to represent the growth of structure in a different cosmol- 
ogy. For example, how to predict the evolution of structure 
in a cosmology with the parameters preferred by analyses 
of the 3rd year WMAP data using the Millennium Simula- 
tion, which was carried out using parameters based on an 
analysis of the 1st year WMAP data. This procedure can 
be repeated for an arbitrary number of target cosmologies, 
fully covering the allowed region of the cosmological param- 
eter space for the standard ACDM model. In §2.1 we provide 
the background and motivation behind our algorithm, and 
in §2.2 we outline a practical implementation of these ideas. 

2.1 Theoretical motivation 

Let us denote by P(k) the linear matter power spectrum ex- 
trapolated to z = in our original cosmology (WMAP1 in 
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the example we give below). The linear theory power spec- 
trum at redshift z is then D 2 (z)P(k), where D(z) is the 
linear growth factor for this original cosmology in units of 
its present-day value, i.e. D(0) = 1. We can now define the 
variance of the linear density field as a function of smoothing 
scale through: 



a (R, z) = D" 



_l r 

4Wo 



dkk 2 P(k)W 2 (kR) = D 2 a 2 (R) (1) 



where R is a comoving smoothing scale, and W(kR) = 
3(sin(fcJ?) — kR cos(kR)) / (kR) s is the Fourier transform of 
the real-space top-hat filter. 

For a specific realisation of P(k), the approximation 
of lZel'dovichl (|l97Cj ) relates Eulerian comoving positions x 
at redshift z to the corresponding initial, or Lagrangian, 
comoving positions q through: 



x(q,z) =q- D(z)S{q), 



(2) 



and in this approximation, the peculiar velocities at redshift 
z are given by: 



v(q,z) 



l + z 



S(q) 



(3) 



The dots here denote derivatives with respect to time and 
S(q) is the linear displacement field of the realisation ex- 
trapolated to z = 0. 

Let us now consider a "target" cosmological model 
(WMAP3 in our example below) which we seek to study 
using simulation results obtained within our original cos- 
mological model. In the following, unprimed quantities will 
continue to refer to the original cosmology, while primed 
quantities P', D' , a' , a' , S' will refer to the target cosmol- 
ogy. Let us assume Zf to be the latest time at which the 
target cosmology must be matched (often Zt = 0). The first 
step in our procedure is to decide a length scaling from orig- 
inal simulation box-size to target box-size L' — s L. This 
is defined so that the linear fluctuation amplitude a'(R,zt) 
in the target cosmology over the range [Ri , R2] is as close 
as possible to that in the original cosmology over the range 
[s _1 Ri, s _1 R2] at a redshift z*(zi). Thus we minimize: 



^rms (^) 2*) 



1 



\n{R 2 /Ri 
R2 dR 
... ~R 



ais' 1 R,z,) 



<T>(R,zt) 111 

over s and z, (see lZheng et al.|[2002l : lHarker et al.ll2007l . for 
a similar approach). 

In Press-Schechter theory the halo mass function is de- 
termined completely by the variation with scale of <j(R, z) 
the smoothed linear variance of the field. EqfJ] is thus 
equivalent to minimizing the difference between the halo 
mass functions in the target and scaled original cosmolo- 
gies over the (target) mass range [M(i?2), M(Ri)]. This sets 
the range of scales included in Eq(4] M(s~ 1 R2) is the mass 
of the largest halo in the original simulation at z = z t and 
M(s~ 1 Ri) that of the smallest resolved halo. For the top-hat 
filter we are using we have M(R) = ^-po R 3 . 

We take the original simulation at z» after scaling all 



lengths by s, to be our first approximation to a simulation 
of the target cosmology at zt. Earlier redshifts in the target 
cosmology (z > Zf) are then matched to those in the original 
cosmology (z > z*) through: 



D'(z' 



D(z t ) 



D'(z ( ) 



(5) 



If needed, peculiar velocities in the target cosmology 
can be obtained approximately from those in the original 
cosmology through: 



D'(z') {l + z) 

S—r^ L -T 7B, 

D(z) (l + z') 



(6) 



but in practice we use a more accurate procedure which we 
now describe. 

At each redshift in the target cosmology we can define 
Rnl(z') through a'(R nU z') = 1 or a' (R nl ) = 1/D'(z'). The 
fluctuations in the target cosmology on all scales larger than 
R n i (i.e. k < R^ 1 ) are in the linear regime and thus the cor- 
responding displacements and velocities can be corrected for 
residual differences in power spectrum between the two cos- 
mologies. This is done by subtracting the long wavelength 
components of the particle position and velocity fields us- 
ing the Ze'dovich approximation (Eqs. [2] and [3} with the 
following low-pass- filtered displacement field, 



n , ls / -i p v jD(z)S k (q) for \k\ < s/R nl 
D(z)S k {q;s Rm) = < 

10 for |fc| > s/Rni 



(7) 



and then adding them back in after scaling each Fourier 
component of the displacement field by the square root of the 
ratio of the amplitude of target to original power spectrum, 




1 1/2 



S sh (q) for |fc| < V-Rni, 
for |fe| > l/R nl , 



(8) 



which is equivalent to the long-range displacement field in 
the target cosmology. Note that the corresponding veloc- 
ity fields are derived from these displacement fields through 
Eq.|3] Also note that, in Fourier space, the low-pass-filtered 
displacement and velocity fields of the original simulation at 
z can be calculated from those used to generate the initial 
conditions. 

A similar technique has been previously studied 
in the literature, alth o ugh f or a d iffere nt pur pose, by 
iTormen fe Bertschingei] (l996l ) and ICold i|l997l) . These 
authors used linear perturbation theory to add long- 
wavelength modes to a n N-bo dy simulation of limited box- 
size. As noted by IColel l|l997D . a naive implementation ne- 
glects the coupling of these long-wavelength modes to the 
smaller-scale perturbations already included in the simula- 
tion. This produces systematic and large biases in the inter- 
nal properties and clustering of haloes, especially massive 
ones. We note that in our case this problem is considerably 
alleviated since we impose relatively small changes on the 
amplitudes of long-wavelength modes. The correct coupling 
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between large and small scales is included automatically to 
lowest order. 



lemen- 



2.2 Implementation 

We now describe briefly the steps of a practical impler 
tation of the ideas presented in the previous subsection 

1) Calculate the scale factor, s, and the redshift z, by min- 
imizing Eq. [3] 

2) For each output z > z* of the original simulation calcu- 
late the corresponding z (z) in the target cosmology using 
Eq.[U 

Then for each z — z pair: 

3) Calculate 7? n i . 

4) Calculate S(q; s~ R n \) for the original simulation. 

5) Invert x — q — D(z)S(q; s~ 1 R n {) in the original simula- 
tion to obtain q(x\ s~ 1 R TL \, z). This can be done recursively: 



q i =x- D(z)S(q; s 1 R n 



(9) 



with q = x as the starting condition and requiring that 
\q i — q i _ 1 \ < 0.01/i _1 Mpc, for example, for convergence. 

6) Remove low-pass-filtered displacement and velocity 
fields from the original simulation data at redshift z. 



7) 
7.1) 
7.2) 
7.3) 



Scale to the target cosmology: 
Box size: L — > L' = sL. 
Redshift: z — ► z'. 
Particle mass: 



Q'H' 2 L' [ 



7.4) Residual nonlinear velocities: 



(1 + *') 



' 1 1/2 



(l + z)m p L' 
n' m H' 2 L'\l + z') 



1 1/2 



(10) 

(11) 

(12) 



where unprimed and primed cosmological parameters are 
present-day values in the original and target cosmologies, 
respectively. 

8) Compute new long-range displacement and velocity 
fields by scaling the Fourier components of the original fields 
according to Eq. [8] using the power spectra of the two cos- 
mologies. Using the Zel'dovich approximation, add these 
low-pass filtered fields to the scaled positions and velocities 
of the particles. 



3 NUMERICAL SIMULATIONS 

In this section we provide details of the N-body simulations 
used to illustrate and probe the accuracy of our algorithm. 
First, we describe the numerical and cosmological parame- 
ters of our simulations, which were carried out at two dif- 
ferent resolutions and in two different cosmologies. Next, we 
discuss how we apply our algorithm to these simulations to 
test our scaling procedures. 
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Wl 


3.34 x 10 11 


865.16 


0.250 


0.045 


0.73 


0.9 


hWl 


3.34 x 10 s 


86.516 


0.250 


0.045 


0.73 


0.9 


W3 


4.92 x 10 11 


1000.0 


0.238 


0.0416 


0.732 


0.761 


hW3 


4.92 x 10 s 


100.0 


0.238 


0.0416 


0.732 


0.761 



Table 1. The values of some of the basic parameters used in 
the simulations. The columns are as follows: (1) The name of the 
simulation. (2) The mass of a dark matter particle. (3) the side 
of the computational box (4) The matter density. (5) The baryon 
density. (6) The Hubble parameter in units of 100 kms/s/Mpc. 
(7) The linear fluctuation amplitude at z = 0, erg. In both cases, 
the primordial spectral index is n s = 1 and the dark energy is 
assumed to be a cosmological constant. 



3.1 Direct N-body Simulations 

In this paper we compare results from cosmological simula- 
tions at two very different resolutions: large- volume simula- 
tions are used to study the performance of our algorithm in 
reproducing large-scale structure; smaller volume but higher 
resolution simulations are used to study its performance in 
reproducing the properties of individual dark matter haloes. 
The latter simulations are also used as input to a semi- 
analytic model for galaxy formation which allow us to com- 
pare the properties of the galaxies predicted in our scaled 
and true simulations. 

The first of our two cosmologies is that of the Millen- 
nium Simulation l|Springel et al.ll2005l ) and is consistent with 
an analysis of the 1st y ear of data from the WMAP satel- 
lite |Spergel et al.|[2003t ). The corresponding large and small 
volume simulation are labelled as Wl and hWl, respectively. 
The second cosmol ogy corresponds to an analysis of the 3r d 
year WMAP data l|Spergel et al.ll2007i : [sirichez et al.ll2006T l. 
and we refer to the large and small volume simulations as W3 
and hW3, respectively. The specific values of the parameters 
for the two models can be found in Table 1. 

Each of our four simulations contains 512 dark matter 
particles. For the W3 simulation the size of the box is L' = 
1000/i _1 Mpc and the redshifts at which we store the particle 
data are z' = (0,0.5, 1.0,2.0). For Wl the corresponding val- 
ues are L = 865.164/i _1 Mpc and z = (0.57, 1.09, 1.68, 2.92). 
These figures were set so that at each output the respec- 
tive box sizes and linear amplitude of fluctuations minimize 
Eq.Has described in §2.1 for [Ri,R 2 ] = [1, 10]/l _1 Mpc. The 
hWl and hW3 simulations have the same cosmological param- 
eters as Wl and W3, but have a 10 times smaller box and 1000 
times smaller particle masses. 

The starting redshifts are z = 60.0 and z — 78.08 
for the Wl and W3 simulations, respectively. The Plummer- 
equivalent softening length has the same value in units of the 
respective box sizes, e = L/10000. The linear power spec- 
trum is calculated for each cosmology using the CAMB package 
l|Lewis et al.ll2000h . The initial density fields are generated 
using the Zel'dovich approximation to perturb a set of par- 
ticles initially distributed on a cubic grid. Finally, all simu- 
lations were e volved with a m emory efficient version of the 
Gadget2 code (|Springelll20"05l ). 

We note that the random number sequence used to set 
initial perturbation amplitudes and phases was the same in 
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Figure 1. The linear mass variance cr 2 (M) smoothed with a real- 
space top-hat filter as a function of the smoothing radius for var- 
ious cosmologies. Symbols in the upper panel show the expected 
values for our standard WMAP1 (black squares) and WMAP3 
(red triangles) cosmologies, and for two cosmologies where all 
WMAP3 parameters are retained except that Q m is raised to 0.4 
(blue triangles) and reduced to 0.1 (green triangles). The con- 
tinuous lines all derive from the WMAP1 cosmology but in the 
red, green and blue cases our scaling in length and our reassign- 
ments of time have been applied in order to match the corre- 
sponding "WMAP3" model as well as possible over the range 
lh _1 Mpc < R < 10/i~ 1 Mpc. The quality of these fits is bet- 
ter appreciated in the bottom panel which shows the difference 
between symbols and curves as a percentage. 

the two different cosmologies simulated at each resolution. 
This allows us to study subtle differences in the final halo 
and matter distributions and provides stringent and clean 
tests for our algorithm. 

3.2 Scaled Simulations: the algorithm in action 

In this section we give details of how the hWl and Wl simula- 
tions were scaled and perturbed, according to the algorithm 
described in §2, to predict the dark matter clustering and 
the halo properties expected in the cosmological model of 
the W3 and hW3 simulations. 

The first step in our procedure is to scale comoving po- 
sitions and to re-label redshifts. As discussed in §2.1, this 
requires finding a new box-size L' and redshift z, for the 
outputs of hWl and Wl that best reproduce the linear mass 
variance at z — of the target WMAP3 cosmology on re- 
solved nonlinear mass scales. 

Linear variance is displayed as a function of scale for 
the original (WMAP1) and target (WMAP3) cosmologies 
in the upper panel of Fig. [T] together with the best fit of 
the WMAP1 variance to the WMAP3 variance after scaling 
lengths and times (see Eq. [4} . The range of smoothing radii 
used in the fit is [1 - 10]/i _1 Mpc in the WMAP3 cosmology, 
corresponding to [10 11 - 10 15 ] hT 1 M for our assumed top- 



Figure 2. Constraints on s and D(z t ) obtained by minimiz- 
ing the logarithmic difference between the WMAP3 and scaled 
WMAP1 linear variances over the scale range corresponding to 
resolved haloes in the WMAP3 cosmology. Contours displayed by 
solid, dotted and dashed lines indicate the regions in parameter 
space that produce rms deviations smaller than 1, 10 and 20%. 
over this mass rangee (see bottom panel of Fig. IT) . 

hat filter in real space. This is the halo mass range expected 
to host the galaxies observable in future spectroscopic sur- 
veys. 

In the bottom panel of Fig. [T] we can see the quality 
of the fit. Over the scales included in the minimization, the 
mass variance in the scaled cosmology differs from that in 
the target cosmology by less than 1%. This deviation implies 
a less than 5% difference in the expected number of "3er" 
haloes between the target simulation and the scaled one, 
assumi ng that the mass function of haloes is well described 
by the iPress fc Schechterl (|l974h model. 

The potential of our algorithm is further illustrated by 
the blue and red lines in Fig. [T] which show the results of 
scaling the WMAP1 variance to match the variance in two 
more distant cosmologies. Q m is set to 0.1 and 0.4 in the red 
and blue cases, respectively, with the remaining parameters 
kept at their WMAP3 values. For both cases our scaling of 
WMAP1 reproduces the linear variance of the "WMAP3" 
model satisfactorily and at a level comparable to the fit of 
the real WMAP3 cosmology, although, the scaling of box 
size and the relabeling of redshift is more drastic, as ex- 
pected. The values for s are 3.47 and 0.536, and for D(z,) 
are 0.416 and 1.311 for lower and higher f2 m , respectively. 

The contour lines in Fig.[2]show the constraints on s and 
on D(z t ) for our standard WMAP1-WMAP3 scaling. The 
best fit values are s — 1.156 and D(z t ) = 0.764 (correspond- 
ing to z* = 0.57). The solid, dotted and dashed lines enclose 
regions in this two-dimensional parameter space which give 
rms differences below 1, 10, and 20% respectively, between 
the scaled WMAP1 variance and the WMAP3 variance over 
the range \R\, R2] (see bottom panel of Fig. [1]). 
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Fig. [2] illustrates the degeneracy between these two pa- 
rameters. The banana-like shape of the contours is, in fact, 
useful because it gives some latitude in choosing which s*, 
in the original cosmology should correspond to z — in the 
target cosmology, allowing one to choose z* to correspond 
to one of the stored outputs. 

For the second part of our algorithm we have pro- 
duced low-pass-filtered displacement and velocity fields for 
the WMAP1 and WMAP3 cosmologies on a 512 3 grid a.s 
described in §2.2. As discussed previously, at each redshift 
we included only the contribution of linear modes with 
k < R~l- For z = this criterion implies k < 0.23 ZiMpc -1 , 
while at z = 2 it implies k < 1.01 ftlvlpc -1 . This corre- 
sponds approximately to smoothing the displacement fields 
with top-hat filters of radii 4.3/i _1 Mpc and 0.98 /i _1 Mpc, 
respectively. 

For the specific simulations we consider in this paper, 
our results imply that if we expand the z = 0.57 outputs of 
Wl and hWl by a factor of 1.156, and then apply the com- 
puted large-scale corrections, the results should be almost 
identical to the W3 and hW3 simulations at z = 0. Below 
we refer to these resulting scaled catalogues as W3' and hW3' 
and we examine in detail how close they are to W3 and hW3, 
respectively. 



4 COMPARISON WITH DIRECT N-BODY 
SIMULATIONS 

In this section we assess the accuracy of our procedure by 
comparing the properties of dark matter, halo and galaxy 
catalogues derived from the scaled W3' and hW3' simulations 
with similar catalogues derived from the true W3 and hW3 
simulations. 



4.1 Mass clustering 

We start testing our algorithm by examining the spatial 
distribution of dark matter particles. Fig. [3] shows the di- 
mensionless power spectra, A 2 = k 3 P(k) /2ty 2 , in real and 
redshift space, measured in the W3 simulation and in the Wl 
simulation after three different scaling procedures have been 
applied. (The power spectra were computed by assigning 
particles onto a 512 3 grid using the Cloud-in-Cell scheme, 
and then Fast Fourier Transforming.) 

The first scaling (red lines) consists of adjusting, in a 
scale-independent way, the linear amplitude of fluctuations 
in spheres of 8/i -1 Mpc of the z = 0.57 output of Wl to match 
the fluctuations of the z = output of W3. The second scaling 
(green lines) includes a rescaling of the box size in addition 
to the relabeling of the time variable as described in §3.2. 
In these two cases, the velocities were scaled according to 
Eq. [T5] The third scaling (blue lines) uses the linear-theory 
corrections on large scales and thus represents our algorithm 
in full. 

As can be appreciated in the top panel of Fig. [3] all the 
described procedures perform moderately well, at least on 
real space. However, the closer look given in the lower sub- 
panels reveals the deficiencies of each of them. For the sim- 
ple amplitude scaling, panel (b) shows that differences reach 
20% even on scales comparable to the fundamental mode of 
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Figure 3. Comparison of the dimensionless power spectrum, 
A 2 (k) = k s P(k) /2ir 2 , measured in the W3 simulation with those 
of the Wl simulation after different types of scaling. Note that in 
panel (a) the redshift-space spectra have been offset by 1 dex for 
clarity. Poisson shot-noise level. Black and brown lines show the 
power spectra of the dark matter in the W3 and Wl simulations, 
respectively (note that we show the Wl simulation at Z = 0.57 and 
that the black dashed line lies behind the dashed blue line). Red 
lines show the result of scaling up the amplitude of the Wl simu- 
lation at 2 = 0.57 to match the linear amplitude of perturbations 
in spheres of 8h~ x Mpc. Green lines display the spectrum after 
scaling both the simulation box and the redshift. In both cases, 
velocities have been scaled according to Eq. [6] Finally, blue lines 
illustrate the result of including our linear theory correction in 
addition to the length and redshift scalings. The deviations af- 
ter applying the different procedures can be better appreciated 
in panels (b), (c) and (d), which show the fractional difference 
^wi W/Aj^fc) — 1. In each of the four panels, power spectra com- 
puted in real and redshift space are displayed as solid and dashed 
lines respectively. 



the simulation box. As the phases in the two simulations co- 
incide in units of box size (and not in units of ft _1 Mpc), the 
ratio is very noisy. As shown in panel (c), these stochastic 
fluctuations disappear when our length and time scalings are 
applied. Nevertheless, a strong disagreement persists and an 
overestimation of nearly 10% remains in the real-space and 
redshift-space clustering on large scales. These differences 
result from the different overall shape of the linear theory 
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Figure 4. The distribution of the differences in x-coordinate for 
dark matter particles before (dashed lines) and after (solid lines) 
our large-scale correction is applied. Different colours show results 
at different redshifts. 



spectra of the Wl and W3 simulations, and from the fact that 
we have normalised them to agree on nonlinear scales. In 
addition, the residuals in panel (c) show an oscillatory pat- 
tern. This is a consequence of the different BAO wavelength 
in the target and scaled simulations, a reflection in part of 
the different fib, but mostly of the rescaling in length. The 
last panel (d) demonstrates that these large-scale features 
are efficiently removed by our linear correction procedure. 

As can be seen also in panels (b) and (c) of Fig. [3] 
the velocity scaling of Eq. 6 does an extremely good job 
in matching streaming motions (and thus the redshift-space 
power spectra) on large scales, but fails for k > 0.1/iMpc - 
because it is inappropriate for nonlinear quasi- virialised mo- 
tions. Panel (d) shows that our preferred correction on these 
scales (Eq. 12) does a much better job. 

As we show in Fig. [4] after our large-scale correction 
the ID rms difference between the position of particles in 
the true and "scaled" simulations is 0.038/i _1 Mpc at z — 0. 
This is a factor of 5 smaller that the rms difference if only 
the length and time scalings are performed (0.197/i _1 Mpc). 
At higher redshifts the improvement provided by our large- 
scale correction is even greater. For instance, at z = 2 it 
reduces the differences in each dimension by a factor of 10; 
the rms values are 0.072 /i _1 Mpc and 0.0073 /i _1 Mpc before 
and after correction, respectively. 

We examine the performance of our full scaling scheme 
at different redshifts in Fig. [5] Over the range of scales dis- 
played, the dark matter power spectrum of the W3 and W3' 
simulations are in very good agreement. At almost all red- 
shifts the deviations on very large scales (A: < 0.1/iMpc _1 ) 
are of order 0.2% or less. On smaller scales the amplitude 
is systematically lower in the scaled dark matter distribu- 
tion, a discrepancy which gets larger at lower redshift. This 
small difference (5% in real space and 10% in redshift space 
at k ~ lftMpc -1 ) is due to differences in the halo mass- 
concentration relations. In the standard picture for structure 
formation, the concentration of a halo is related to the mean 
density of the Universe at the time of its formation (which 
is usually defined in terms of its mass accretion history). 
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Figure 5. Relative difference, A^ 3 , (fc)/Ay 3 (fc) — 1, between the 
dimensionless dark matter power spectra measured in the scaled 
and true simulations as a function of redshift. Solid lines com- 
pare W3' and W3 and dashed lines hW3' and hW3. Black curves 
refer to real-space and red to redshift-space power spectrum. In 
each panel, the arrow indicates the wavelength, k n i, that satisfies 
A 2 (fc n ;) = 1. This is the largest fc-value included in our linear the- 
ory corrections. The legend in each panel indicates the redshifts 
of the original (Wl) and target (W3) simulation outputs. 



Even though our true and scaled cosmologies have, by con- 
struction, the same perturbation amplitudes at matched red- 
shifts, their growth histories differ when expressed in terms 
of redshift. Thus the mean collapse redshift of haloes of given 
mass will also be different l|Dolag et alj|2004l ; Ijennings et al.l 
2009). As we will see below, this implies that haloes in the 
scaled simulations are less concentrated than in the true 
simulations. The effect increases with decreasing redshift 
reaching 10% by z — 0. On the other hand, the absence 
of wiggles in the residual around the BAO scale implies that 
our Zel'dovich-based corrections do appropriately compen- 
sate not only for the linear signal but also, at lowest order, 
for its nonlinear degradation. 

A visual impression of the quality of the scaling is given 
by Fig. [6] This figure shows the dark matter density field in 
real space at z = measured in 50/i -1 Mpc slices through the 
W3 (left column) and W3' (central column) simulations. Un- 
derdense regions are lighter while clusters and filaments are 
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Figure 6. The mass density field in a slice of thickness 50h~ x Mpc and side 1000/i —1 Mpc at z = for the W3 (left column) and the 
W3' simulations (middle panel). The difference between these two maps is shown at greatly enhanced contrast in the right column. From 
top to bottom, the dark matter field has been smoothed with Gaussians of radius 2, 5 and 20/i -1 Mpc. The intensity of each pixel is 
proportional of the logarithm of the density. 
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Figure 7. Comparison of the comoving number density of FoF 
haloes. We plot the multiplicity function p _1 M 2 dA r /dAf where 
p is the mean density of the Universe. The curves, from top to 
bottom, represent redshifts z = 0, 0.5 and 1.0. Solid lines dis- 
play measurements from the direct simulations hW3 and W3 while 
dashed lines show results from the scaled catalogues W3' and hW3' . 
Vertical dotted lines indicate the 100 particles limit for the hW3 
and W3 simulations. The lower panel shows the difference "scaled" - 
true as a percentage. 



darker The differences between the two maps are displayed 
with enhanced contrast in the right column. The density 
fields have been smoothed with Gaussian kernels of radius 
ro = 2, 5 and 20 ft. - Mpc (top, middle and bottom rows 
respectively). In this way we accentuate the differences be- 
tween the dark matter fields on different scales. 

In every row, we can see that the two simulations have 
identical structures - voids, filaments and clusters are found 
at exactly the same locations. Small scales, best seen in the 
upper panel, show a slight positive bias in dense regions 
which results from the higher concentration of dark mat- 
ter haloes in the true W3 simulation than in the scaled W3' 
simulation. Small shifts in the position of clusters appear as 
black-white dipoles in the right column. The discrepancies 
between the maps are small in all cases and are strongly re- 
duced when we smooth more strongly. In the bottom row the 
two fields are virtually indistinguishable. This is consistent 
with our power spectrum analysis: the matching procedure 
performs remarkably well on large scales but does not per- 
fectly reproduce the small-scale nonlinear clustering. 



4.2 The properties of DM haloes 

We now compare the distribution of dark matter haloes in 
the true (W3) and scaled W3' simula tions. These catal ogues 
were built using a FoF group finder (|Davis et alJll985h with 
a linking length equal to 20% of the mean interparticle sep- 
aration, which approximately selects objects with 200 times 
the mean density. At z — the number of haloes with at 



least 50 particles in the W3 and W3' catalogues is 138, 284 and 
131, 356, respectively. These figures are 109, 651 and 109, 176 
for the hW3 and hW3' catalogues. 

Fig. [7] shows the abundance of dark matter haloes as 
a function of mass in our two simulations. Over the mass 
range covered by the catalogues, the mass functions are in 
very good agreement. Differences are below 10% over 4 or- 
ders of magnitude in mass and over 7 orders of magnitude 
in number density! However, focusing on the bottom panel, 
where differences between the curves can be seen in more 
detail, it is evident that there are significant systematic dis- 
crepancies. 

The scaled simulations W3' and hW3' seem to underpre- 
dict the abundance of haloes compared to the W3 and hW3 
simulations, or equivalently, FoF masses are systematically 
higher in the true simulations than in the scaled simulations. 
In particular, the difference in abundance seems larger for 
high-mass haloes, this explains the fact that the total num- 
ber of haloes differs by 5% between W3 and W3' , but only by 
0.5% between hW3 and hW3'. There are a number of possible 
reasons for this disagreement. A ~ 5% disagreement is ex- 
pected due to the residuals in the mass variance fitting (See 
Fig- HI. In addition, a fixed linking length identifies objects 
with lower mean overdensities relative to the mean density 
if their concentrations are lower, therefore underestimating 
the mass of objects identified at low redshift. Finally, an im- 
portant assumption of our work is that the halo mass func- 
tion is universal (i.e. independent of redshift, cosmological 
parameters, etc) when the masses are expressed in terms of 
5 c /a(M), as predicted by Press-Schechter theory. However, 
recent simulation results suggest that this picture may work 
no better than 10%, and that properties of the density field 
other than the shape of the linear vari ance as a function 
of smoothing radiu s may be important (|Tinker et al.ll2008l ; 
iManera et~aUl2009F) . 

In order to provide a more detailed comparison between 
our catalogues we have matched individual haloes in W3' and 
W3. We do this by selecting all the DM particles that belong 
to each halo with more than 50 particles in W3 and then 
finding the halo that contains most of these particles in W3'. 
Due to statistical fluctuations not every halo in the W3 (hW3) 
catalogue can be matched to one in the W3' (hW3') cata- 
logue. Nevertheless, of the 138, 284 (109, 651) haloes with 
more than 50 particles in W3 (hW3), only 16 (21) are un- 
matched to a halo with at least 20 particles in W3' (hW3'). 

The object-by-object comparison that results from this 
matching is presented in Fig. [5] In this plot we present the 
distribution of differences in the position, velocity, mass and 
concentration of matched FoF haloes for the W3-W3' (dotted 
lines) and for the hW3-hW3' catalogues (solid lines). We show 
the results at z = and for haloes in 6 different mass ranges, 
but we have checked that our results behave in a similar way 
at different redshifts, and for other mass ranges. 

In the top left panel we can see that the centre of 
mass positions of the FoF haloes, agree remarkably well. 
For each of the 3 spatial coordinates, discrepancies are typ- 
ically 85 ft _1 kpc between the l/i -1 Gpc simulation, but only 
30 /i -1 kpc between the 100/i~ 1 Mpc simulations. These un- 
certainties are less than half the typical virial radii of the 
haloes. This agrees with the small differences seen in Fig|5] 
We note that these differences increase roughly by a fac- 
tor of two if our large-scale correction is not included. In 
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Figure 8. The distribution of the differences in position, velocity modulus, mass and concentration for FoF haloes at z = in W3 
and W3' (dotted lines) and in hW3 and hW3' (solid lines). Results are shown for 6 disjoint mass ranges. Black, red and green curves show 
distributions for haloes with [200,700], [700,1500] and > 1500 particles in W3 and hW3. The standard deviation of the distributions is 
given in the legend of each panel. 



panel (b) we display differences in the ID velocity modulus 
(V 2 — V 2 + Vy + V 2 ) of the centre of mass of our haloes. 
We find that this quantity is also matched very accurately: 
better than ~ 5% rms, and with no measurable bias for both 
pairs of simulations in all the mass bins inspected. 

On the other hand, the FoF masses are slightly less 
precise. The uncertainty in our scaled halo mass is just under 
10% with a systematic underestimation of ~ 5%. We also 
note that the accuracy at which we match this is greater for 
our high-resolution simulations. 

In contrast, in panel (d), we compare concentrations of 



FoF haloes using the ratio C = Vmax/Vzoo, where Vmax is 
the maximum circular velocity of the main subhalo of the 
FoF halo and V200 is the circular velocity at -R200 , the radius 
enclosing a mean density 200 times the critical value. The 
discrepancies here are typically larger in the high-resolution 
pair, although still under 10% with a bias of about 5%. The 
lower mean concentration of haloes in W3' is responsible for 
the underestimation of the real-space power spectra on small 
scales seen in Fig. [5] We expect this to impact the spatial 
distribution of galaxies, if at all, only on very small scales. 

In Fig. [9] we show the accuracy with which our scaling 
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Figure 9. Comparison of the concentration (left panel) and of the spin parameter (right panel) for haloes identified in the true 
(unprimed) and scaled (primed) simulations at z = 0. For clarity, we have plotted the result only for haloes containing more than 500 
particles. The red lines indicate the result expected if there were no systematic biases between the halo properties in our two simulations 
(i.e. the one-to-one relation), while the blue lines indicate the straight line that minimises the rms perpendicular residuals. 



procedure is able to reproduce two important properties of 
dark matter haloes. In the left panel we compare concen- 
trations (as defined above) of haloes identified in the scaled 
(C") and true simulations (C). In agreement with Fig. [8] 
there is a small bias in the average value which results from 
differing halo growth with redshift in our two simulations: 
the concentration of a halo is underestimated by 0.04 in the 
scaled simulations. In the right panel we compare the spin 
parameter defined as A = J/(V2M2oo V200-R200) where J 
is the halo angular momentum and M 200 is the total mass 
(both within R200) (|Bullock et al.l l2001). This panel demon- 
strates that, albeit with scatter, spins are well reproduced 
by our velocity and position scaling. 

We now explore whether the differences displayed in 
Fig. [8] and Fig. [9] produce a noticeable effect on halo clus- 
tering by comparing the power spectra of the matched cata- 
logues. Fig.[in]presents results in real space (black lines) and 
in redshift space (red lines) at z — 0, 0.5, 1 and 2. As found 
above for the dark matter catalogues, on large scales the 
two halo power spectra are almost identical. In real space 
the differences are less than 0.1% while in redshift space 
they are still below the percent level. This behaviour repeats 
at all redshifts. On smaller scales, haloes in W3' tend to be 
more clustered than their counterparts in W3. Overall, the 
distributions agree, over all relevant scales, within 1%, an 
uncertainty much smaller than those likely to be introduced 
by our poor knowledge of the galaxy formation process. 

It is important to note that in these measurements we 
have used haloes above a fixed mass. As a result, the sam- 
ples displayed at different redshifts correspond to peaks of 
different height in the density field, as illustrated by the bias 
values given in each panel of Fig. 10. Thus, the fact that we 
recover the expected clustering at all redshifts suggests that 



the bias-mass relationship is very well reproduced by our 
algorithm, even for very rare peaks. 



4.3 Models for galaxy formation 

The final step in the creation of galaxy mock catalogues 
is to connect the predictions for dark matter haloes with 
the observable properties of the galaxies which they host. 
For this, a crucial requirement is the construction of accu- 
rate merger histories. Thus we now e xplore how well our 
algorithm performs in this respect (see lNeistein et aill2009l . 
for an attempt to describe merger histories in a cosmology- 
independent way). We do this by comparing the galaxies 
built using hW3' with those built using an identical proce- 
dure on hW3. 

There are several different approaches to create realis- 
tic galaxy populations in the dark matter haloes of N-body 
simulations. Here we use so-called semi-analytic models of 
galaxy formation in which the key physical processes that 
determine the properties of galaxies are followed from very 
early times in a self-consistent way within the hierarchically 
growing halo population. Such semi-analytic modelling rep- 
resents a hard test for our algorithm since the final proper- 
ties of galaxies depend not only on the present-day mass of 
their host haloes but also on their whole assembly histories 
- merger rates, progenitor mass function, formation times, 
spins - and indirectly on the properties of the environment 
that affect these quantities. 

The specific model we use in this pa per is that developed 
by th e Munich group as describ ed in De Lucia fc Blaizotl 
i|2007l ). with the parameters of IWang et al.l (|200Sl ). This 
model roughly reproduces the observed z = K-band lu- 
minosity function (among other galaxy properties) in the 
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Figure 10. Relative difference, A* 3 , (k)/A% 3 (k) - 
dimensionless halo power spectra measured in the scaled and true 
simulations as a function of redshift. We have used only matched 
FoF haloes with at least 50 particles in W3. This limit corresponds 
to a mass of 2.6x 10 12 /i _1 A/ Q in the WMAP3 cosmology. The bias 

of the sample, defined as the mean value of b = A 2 (fe)/A^ m (fc) 
on scales k < 0.1/lMpc — 1 , is shown in the legend of each subpanel 
and is seen to vary strongly with redshift. 
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Figure 11. The K-band luminosity function at z = constructed 
from the hW3 (green squares) and hW3' simulations (red triangles). 
For comparison, we have also included the luminosity function of 
the 2dFGRS galaxies as computed from Cole et al.(2001). 



ist: there is a small underestimate of the number of galaxies 
brighter than a given magnitude in the hW3' catalogue. This 
difference is weakly magnitude-dependent and reaches 5% at 
the faintest magnitudes. It is consistent with the differences 
in the halo mass function seen in Fig. [7] 

A further comparison of the results of the semi-analytic 
modelling is presented in Fig. 1121 In this plot we show the 
distribution of K-band magnitude differences, K hV3 — K hv3 i , 
for the central galaxies of "matched" haloes (c.f. §4.2). The 
predicted properties of galaxies in the WMAP3 cosmology 
are quite accurately reproduced by our scaled simulations. 
Typical random differences are below 0.25 magnitudes over 
the whole magnitude range and systematic errors are below 
0.13 magnitudes. 



WMAP3 cosmology, and is adequate for our tests. For both 
hW3 and hW3' we have generated halo merger trees with the 
same time resolution as in the Millennium Simulation, i.e. 
roughly 300 Myr for z < 1. All levels in the two merger trees 
are set to correspond appropriately after scaling lengths and 
masses and reassigning redshifts. This involves matching 
outputs in the hW3 and the original hWl simulations as above. 

The present-day K-band luminosity functions computed 
by this model are shown in Fig llll The green triangles refer 
to galaxies in the true model hW3, while the blue squares 
refer to galaxies in the scaled model hW3'. For comparison 
we have also inclu ded the observed luminosity function from 
ICole et all (|200lT ) . In both cases our galaxy samples are com- 
plete down to Mk — 5 log h = —18. Above this magnitude 
threshold the catalogues contain 76562 and 71652 galaxies 
for hW3 and hW3', respectively. Overall the agreement is quite 
impressive, especially if one remembers that we are compar- 
ing results for simulations run with significantly different 
sets of cosmological parameters. However, discrepancies ex- 



5 CONCLUSIONS 

In this paper we have introduced an algorithm that allows 
us to scale the output of a cosmological N-body simulation 
carried out for one specific set of cosmological parameters 
(for example the WMAP1 parameter set of the Millennium 
Simulation) so that it faithfully represents the growth of 
structure in a different cosmology (for example, one with 
the lower fluctuation amplitude preferred by analyses of the 
WMAP3 data). Our method builds on the fact that the lin- 
ear fluctuation amplitude as a function of smoothing scale 
is sufficient to determine not only the mass function and 
assembly history of haloes (according to extended Press- 
Schechter theory) but also the ensemble properties of the 
large scale structure (in the Zel'dovich approximation). We 
therefore re-scale the box-size of the simulation and reassign 
redshifts to the outputs so that the shape of the smoothed 
linear power spectrum matches that of the target cosmology 
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Figure 12. The differences in K-band absolute magnitude for 
central galaxies in the hW3 and hW3' catalogues. These galaxies 
represent the ~ 56% of all the galaxies with Mk — 5 log 10 ft < — 18 
in our semi-analytic catalogues. The median of the distribution 
for the sample as a whole is 0.135 while the dispersion is 0.233. 



on the scales which produce resolved dark matter haloes. 
This scaling and redshift reassignment then imply a rescal- 
ing both of particle mass, and of the nonlinear "thermal" 
velocities in each output. Finally, we use the Zel'dovich ap- 
proximation to modify the amplitude of the long-wavelength 
quasilinear modes in order to account for the difference in 
shape on large scales between the smoothed linear power 
spectrum of the original cosmology (after scaling) and that 
of the target cosmology. This modification is effectively a 
rescaling of the displacements and peculiar velocities pro- 
duced by each individual long-wavelength (low k) mode. 

We tested our method by applying it to an N-body sim- 
ulation based on the WMAP1 parameters to mimic the mass 
distribution expected in the WMAP3 cosmology. We com- 
pared the results with those from a matched N-body simula- 
tion run directly in the target cosmology (WMAP3). These 
two cosmological models span the range of power spectrum 
shape and amplitude among the models currently considered 
plausible. Thus we expect our comparison to give conserva- 
tive estimates of the accuracy of our method. 

We found that on large scales the target mass dis- 
tribution is recovered to high precision (< 1% on scales 
k < 0.2 ftMpc -1 ) in both real and redshift space. On smaller 
scales the differences increase ( 5%) and can be attributed 
to differences in the halo mass-concentration relation be- 
tween the true and scaled simulations. These differences are 
smaller when one considers the clustering of dark matter 
haloes, - the distribution actually relevant for galaxy clus- 
tering. In the latter case, differences are at the subpercent 
level over the full range of scales we examined. We also tested 
the algorithm on an object-by-object basis concluding that 
the position of individual FoF haloes is reproduced to bet- 



ter than 90 h 1 kpc, and their velocities and FoF masses to 
better than 10%. 

The dispersion in the matched properties introduces no 
appreciable effects in the overall clustering statistics of the 
halo population. Finally for a representative semi-analytic 
galaxy formation model, the K magnitudes of central galax- 
ies are accurate to 0.25 magnitudes with a systematic error 
of 0.13 magnitudes. 

Understanding the connection between galaxy forma- 
tion and the background cosmological model is an impor- 
tant challenge that cosmologists will face in the near future. 
A key step towards this is the creation of realistic physically 
based theoretical models for the galaxy distribution for a 
wide range of cosmologies. Our method offers a route to this 
goal. Its relatively small computational cost enables explo- 
ration of a multidimensional cosmological parameter space 
using a suitable N-body simulation of a single cosmology. 
In this way the task of creating realistic mock galaxy cata- 
logues for all allowed cosmological parameter sets is greatly 
simplified. It also becomes possible to update existing high- 
resolution simulations to any favoured cosmological model 
at a small fraction of the original computational cost. These 
goals are important not only to understand the physical pro- 
cesses underlying galaxy formation, but also to ensure op- 
timal exploitation of the information encoded in the very 
large (and expensive) survey datasets that are currently be- 
ing assembled. Our scheme makes it possible to study how 
uncertainties in galaxy formation physics will affect our abil- 
ity to constrain cosmological parameters through analysis of 
next-generation galaxy surveys. 
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